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ABSTRACT 

We consider the stability of a configuration consisting of a vertical magnetic 
field in a planar flow on elliptical streamlines in ideal hydromagnetics. In the ab- 
sence of a magnetic field the elliptical flow is universally unstable (the "elliptical 
instability"). We find this universal instability persists in the presence of mag- 
netic fields of arbitrary strength, although the growthrate decreases somewhat. 
We also find further instabilities due to the presence of the magnetic field. One 
of these, a destabilization of Alfven waves, requires the magnetic parameter to 
exceed a certain critical value. A second, involving a mixing of hydrodynamic 
and magnetic modes, occurs for all magnetic-field strengths. These instabilities 
may be important in tidally distorted or otherwise elliptical disks. A disk of finite 
thickness is stable if the magnetic fieldstrength exceeds a critical value, similar 
to the fieldstrength which suppresses the magnetorotational instability. 



Subject headings: accretion discs - instabilities: elliptical, hydromagnetic 
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1. Introduction 

The problem of momentum transport in accretion disks is widely believed to require 
hydrodynamic or hydromagnetic turbulence for its resolution. The origin of this turbulence 
may be sought in the instability of laminar solutions of the equations of hydromagnetics, 
solutions that are compatible with the geometry of accretion disks. The recent history of 
these efforts has taken the form of first recognizing such an instability mechanism, and then 
trying to incorporate that mechanicsm into realistic disk models. 

The magnetorotational instability (MRI) mechanism, originally discovered by Velikhov 
(1959) and Chandrasekhar (1960) and first applied to accretion disks in Balbus & Hawley 
(1991), is of this kind (see Balbus & Hawley (1998) for a review). It appears in rotating, 
magnetized systems in which the specific angular momentum increases outward and in which 
the magnetic field is weak enough that rotational effects are not overwhelmed by magnetic 
tension. 

A second mechanism, that of the elliptical instability considered by Goodman (1993) and 
others (Lubow, Pringle, and Kerswell 1993; Ryu and Goodman 1994; Ryu, Goodman, and 
Vishniac 1996), is also consistent with the accretion-disk setting. This instability mechanism 
has been reviewed by Kerswell (2002). In the setting considered by Goodman et. al., it 
appears to require a secondary in order to enforce departure from rotational symmetry of 
the streamlines via a tidal potential. This is certainly appropriate for binary systems but 
it is likely that, even in the absence of a secondary, the laminar motion in the plane of the 
disk would not be accurately circular, so the elliptical-instability mechanism would appear 
to be a candidate of considerable generality. It does not require a magnetic field. One of 
the conclusions of the present paper is that it further persists in the presence of a magnetic 
field. In the idealized setting of the present problem, the latter may be of arbitrarily large 
strength. However, we also argue that in the setting of a disk geometry, there may indeed 
be a limit on the field strength. 

In this paper we therefore investigate the interaction of a vertical magnetic field with flow 
on elliptical streamlines, on the ground that both magnetic fields and noncircular streamlines 
are likely ingredients in accretion-disk settings. There are similarities with and differences 
from previous work on effect of magnetic fields on the elliptical instability (Kerswell 1994), 
which are discussed in 56. 



- 3- 



2. Formulation 

We consider flow on elliptical streamlines together with a magnetic field and investigate 
linear stability theory. The underlying equations are the Euler equations of fluid dynamics 

u t + u ■ Vu = —Vp + (curl B) x B (1) 

and the induction equation 

Bi + u-VB = BVu. (2) 

We shall assume that div u = and div B = and that the fluid is unbounded. 

It is easy to check that the following steady fields represent a solution of the preceding 
system: 

U = f2(-%,,^si,oY B = (0,0,5), P=%{x\ + x 2 2 ). (3) 

\ 0,2 a l / * 

Here Q and B are constants, and a constant may also be added to the pressure term. More 
general exact solutions of the combined fluid/magnetic equations exist in an unbounded 
domain (Craik 1988); the case in hand is probably the simplest of these. 



2.1. The perturbed system 

Let u, B, p be replaced by U + u, B + b, P + p in equations (1, 2) above, and linearize. 
The resulting perturbation equations are 

u t + U • Vu + u ■ VU = -Vp + (curlb) x B (4) 

and 

b t + U • Vb = B • Vu + b • VU, (5) 

together with the conditions that u and b be solenoidal. These equations allow rotating-wave 
solutions of the form 

u = v (t) exp % (k (t) , x) , b = w (t) exp % (k (t) , x) , p = (t) exp i (k (t) , x) (6) 

where the expression (k, x) denotes the inner product. Because u and b are solenoidal, the 
conditions 

(k(f),v(f)) = 0, (k(f),w(f))=0 (7) 

must be satisfied. 
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Write 

/ —E 

U = Ax where A = ft , | , E = a 1 /a 2 . (8) 

\ 

Then on substituting the rotating-wave expressions from equation (6) into the perturbation 
equations one finds 

k = -A% (9) 
v = -Av - ik(f) + iB (k x w) x e 3 (10) 
w = i(k 3 B)v + Aw. (11) 

Equation (9) can be solved to give 

k = (k cos (ftt -x),Ek sin (Sit - x) , h) , (12) 

where k, ks and x are constants. The pressure coefficient can be eliminated with the aid 
of the solenoidal condition (7). One finds 

-i{<f) + Bw 3 ) = 2k' 2 (A%v) . (13) 

The equation for v now takes the form 

v = C (t) v + % (k 3 B) w, (14) 

where 

/ -E~ x k x k 2 Ek\ \ / E 0\ 

C (t) = -2 (ft/A; 2 ) —E~ l k\ Ek x k 2 +ft —E' 1 0. (15) 
\ -E~ l k 3 k 2 Ek x k 3 / \ / 

It's convenient to break the six-dimensional system consisting of equations (14) and (11) 
into two, one of size four and the other of size two: 

vi = -(2n/k 2 )k 1 (Ek 1 v 2 -E~ 1 k 2 v 1 ) +QEv 2 + imw 1 , 

v 2 = — (2ft//c 2 ) k 2 (Ek±v 2 — E^k^i) — VtE~ l v\ + imw 2 , 

w 1 = imvi — QEw 2 , 

w 2 = imv 2 + VtE^wi, 

where m = k 3 B. These four equations are self-contained and the remaining equations, 

v 3 = - (2ft/£; 2 ) k 3 {Ekiv 2 - E^fovi) + imw 3 , (16) 
w 3 = imv 3 (17) 
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may be integrated once the expression 

ci = Ek x v 2 - E- 1 k 2 v 1 (18) 

is found by solving the four- dimensional system above. Equations (9), (10), (11) and (13) 
imply that 

— (k, v) = im (k, w) and — (k, w) = im (k, v) . 

Thus in solving this system we need to impose the conditions that these inner products are 
zero initially; this will thereafter maintain the incompressibility conditions (7). 

The incompressibility condition provides an alternative way of finding t> 3 and w 3 once 
the equations for v±,v 2 ,wi,w 2 have been solved, provided that k 3 ^ 0. The only cases for 
which k 3 can vanish are those for which the combinations fci^i + k 2 v 2 and kiWi + k 2 w 2 are 
also found to vanish on solving the four- dimensional system above. It is not difficult to 
show that there can be no instability associated with such a solution (see in particular the 
equivalent system (20) below). Accordingly, we henceforth consider only perturbations with 
vertical wave number k 3 ^ 0. 



2.2. Change of variables 

We change to new variables to facilitate subsequent calculations 1 . 

ci = Ek x v 2 - E~ x k 2 vx, 

c 2 = k x vi + k 2 v 2 (= -/C3V3) , 

c 3 = Ekiw 2 - E^k^x, 

c 4 = k x w l + k 2 w 2 (= -k 3 w 3 ) . 



(19) 



This is a time-dependent (periodic) change of variables since k is periodic in t. The equations 
to be solved take the form 

c = D(t)c (20) 
in these variables, with (we have put Q = 1 here to agree with earlier conventions) 



D(t) 



/ 


—2 (E — E~ l ) k~ 2 k\k 2 


-2 


im 





\ 












im 






im 













V 





im 








/ 



(21) 



The origin of this change of variables is related to the existence of the rotating- wave solutions. 
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2.3. General Considerations 

The coefficients of the matrix D depend on the phase angle x appearing in the ex- 
pressions (12) above for the wave vector k. For purposes of studying stability, we may set 
X = 2 . This is easily seen by making the substitution t' = t — x, which eliminates x from 
the equation. For the remainder of this work we take x — 0- 

The system (20) presents a Floquet problem (cf. Yakubovich and Starzhinsii (1975)): 
the stability of the chosen steady solution U, B depends on whether there are solutions of this 
system that grow exponentially with time. This is settled by finding the Floquet multiplier 
matrix M. The latter is defined as follows. Let $ (t) be the fundamental matrix solution of 
equation (20) that reduces to the identity at t — 0. Then, since the periodicity of D is 2tc, 
M = $ (27r) . If any eigenvalue A of M has modulus exceeding one, this implies that there is 
indeed an exponentially growing solution. 

It is familiar in conservative problems that 

Proposition 1 Whenever A is an eigenvalue of the Floquet matrix, so also are its inverse 
A -1 and its complex conjugate A. 

The first statement of this proposition is a typically a consequence of canonical Hamiltonian 
structure, the second a consequence of the reality of the underlying problem. However, the 
system (20) is not canonical, and the matrix appearing in it is not real. We can nevertheless 
establish these familiar properties of the eigenvalues directly from the system (20), as follows. 
The time-reversal invariance of the physical problem is reflected in the existence of a reversing 
symmetry R = diag (1,-1,-1,1) of the matrix D above: RD (— t) = —D (t) R, implying 
that whenever c (t) is a solution so also is Rc (—t). Since the solutions of the Floquet problem 
have the structure c(t) = pit) exp (at), there must also be a solution p(—t) exp (—crt) . But 
the eigenvalues of the Floquet multiplier matrix M are the values A = exp (2na). This shows 
that if A is an eigenvalue of M, so also is exp (— 2na) = A -1 . 

Similarly, under matrix transformation S = diag (1, 1,-1,-1), D goes to its complex 
conjugate. This shows that M = SMS" 1 , i.e., that M and its conjugate have the same 
eigenvalues. 

Immediate consequences of Proposition 1 are the following: first, in the stable case, 
eigenvalues of M lie on the unit circle; second, if, as parameters change, an eigenvalue is 



2 On the other hand, for purposes of solving the initial- value problem, which involves integrating over 
initial wave vectors, we would need to retain it. 
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at the onset of instability, it must have multiplicity two (or higher). The latter conclusion 
is because, in the complex A-plane, the dangerous eigenvalue A must leave the unit circle 
simultaneously with A , which lies along the same ray as A and therefore coincides with 
it when they both lie on the unit circle. Thus a necessary condition for the onset of linear 
instability is a resonance where two Floquet multipliers coincide. 

2.4. Parameters 

There are three dimensionless parameters that figure in this problem. We call them e, //, 
and r\: 

e = ^(E- E~ l ) , /I = k 3 /k , and rj = k B. (22) 

Thus e represents the departure of the streamlines of the unperturbed flow from axial sym- 
metry. In these equations k = a/ ' k 2 + k\ and represents the length of the wave vector if 
e = 0. The magnetic parameter r\ depends not only on the strength of the unperturbed 
magnetic field but also on the wavelength of the perturbation. 

In the matrix D above, the magnetic field enters through the parameter m = k 3 B (which 
is a measure of the magnetic tension force), and we shall continue to use this notation for the 
present, on the understanding that m — /irj. It is also clear that we can use E rather than e 
to measure the departure from rotational symmetry, and we shall do this in some cases. 

3. Analysis 

The Floquet matrix is 

M(e,n,r)) = $ (2tt, e, n, rj) . 

One could map out the stability and instability regions in the e/irj parameter space numer- 
ically by integrating the system (20) systematically for many values of these parameters. 
We in fact do this for a selection of parameter values in §4 below. However, in this and 
the following section we present the outlines and results of an asymptotic analysis based on 
regarding e as a small parameter (details are presented in the Appendices). This is more 
revealing than the numerical reults on their own. It is also of considerable importance in 
interpreting the numerical results, and is quite accurate even for values of e that are not very 
small (cf. Figure 1 below). The calculation proceeds in two steps; finding M, and calculating 
its eigenvalues. 
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3.1. The Floquet Matrix 



The asymptotic analysis is facilitated by the circumstance that, if we put e = 0, the 
coefficient matrix Dq (say) of equation (21) becomes constant. We find the eigenvalues and 
eigenvectors of D in Appendix A.l. There are two complex conjugate pairs of modes. One 
pair reduces to the ordinary hydrodynamic modes in the limit rj — > 0. The second pair are 
magnetic modes with zero frequency at rj = 0. In the weak field limit, the ratio of magnetic 
to kinetic energy is 0(rf/A) for the modified hydrodynamic modes and 0(A/rf ) for the 
magnetic modes. If 77 ^> 1, the kinetic and magnetic energies are near equipartition for both 
types of mode. 

We now turn to a brief description of the asymptotic (or perturbation) procedure. In 
what follows the parameter r\ will be held fixed so, to simplify the notation, we suppress the 
dependence of the Floquet matrix on this parameter: M — M (e, fx) — $ (2ir, e, //). We shall 
need the Taylor expansion 



where the dots indicate higher-order terms in e and fx — /x . The reason for allowing variations 
with fx as well as variations with e is that the region in the e/x plane where instability occurs 
is typically a wedge with apex at a point (e, fx) = (0, fxo) and boundaries fx = (xq + v±e, where 
the slopes v + and z/_ are to be found, (cf. Figure 1 below). We will therefore consider fx of 
the form 



where v may be regarded as a fixed parameter to be chosen later. This will lead us to the 
widest of wedges in the efx plane, of width O (e), excluding other wedges of width O (e m ) 
with m > 2. These higher-order wedges typically occupy a tiny fraction of the parameter 
space (see Figure 4 below). As a result of the representation (24) we may write 



M{e,fx)=M(0,fi )+M e (0,fx )e + M fM (0,fi ) (fx - fx ) + ■ ■ ■ , 



(23) 



fx = fx + ve H 



(24) 



M = M + M : < H 



(25) 



where 



M = M (0, /i ) and M l = M e (0, fx ) + vM^ (0, fx ) . 



(26) 



The matrix Mo is given in Appendix A.l, together with an expression for M e . We effect 
a change of basis, which diagonalizes Mq and simplifies the calculation of M e , in Appendix 
A. 2. The elements of Mi are calculated in Appendix A. 4. 
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3.2. The characteristic polynomial 

Denote by 

p(A,e) = \M(e) - XI\ (27) 

the characteristic polynomial of the Floquet multiplier matrix M (e) and its roots by Ai, A2, • • • 
A necessary condition for stability is that each root lie on the unit circle. To explore the on- 
set of instability for various parameter values, we shall obtain the characteristic polynomial 
p in the form 

p(A,e) =p (A) +Pi(\)e + p 2 (\)e 2 + ••• (28) 

and exploit our knowledge of the roots of po to obtain the roots of p (A, e) in the form of a 
Puiseux expansion in e (Hille 1972)). 

The nature of this expansion depends on the multiplicities of the roots {A^} of po, the 
characteristic polynomial of the unperturbed Floquet matrix Mq. These roots are given by 
the expressions A& = exp27ro"fc where the {a k } are the eigenvalues of the matrix D given in 
Appendix A.l (equation A3) ; they are all distinct. However, it is possible for the multipliers 
{Afc} to be repeated even when, as in the present case, the {a k } are distinct: if — a\ — ik 
for an integer k 7^ 0, then \ k = \i. 

A necessary condition for the onset of instability is that there be a double (or higher) 
root of the characteristic equation, and we henceforth restrict consideration to the case of 
double roots 3 . For definiteness, we suppose A 2 = Ai. Then the Puiseux expansion takes the 
form 

Ax(e) = Ax + e 1 ' 2 (3 1/2 + eft + • • • . (29) 

Substituting into the characteristic equation (and taking into account that p and p' Q both 
vanish at Ai) yields for the coefficient j3i/ 2 of the leading-order correction the equation 

ftp = -2 Pl (Ai) Ipl (AO . (30) 

The two values of /3i/2 give the generic expressions for the change in a double eigenvalue, 
yielding a pair of roots branching from the double root Ai. However, if p 1 (Ai) = 0, this 
expression is inadequate and one must proceed to the next term in order to determine the 
effect of the perturbation on the stability. Under the present assumptions, it is indeed the 
case that p\ vanishes at Ai, as we show in Appendix A.3 4 . 



^Higher order zeros are not ruled out in this problem, since there are three independent parameters, but 
we do not pursue this here. 

4 This "nongcneric" behavior can be traced to the circumstance that the perturbation expansion takes 
place at a codimension-two point, i.e., where two relations must hold among the parameters. 
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We must suppose then that the expansion of p (A, e) is carried out to second order in e: 

p (A, e) = Po (A) + e Pl (A) + e 2 p 2 (A) H (31) 

Then in the Puiseux expansion above /3y 2 — so Ai = Ai + foe + • • • , and f3\ is found by 
solving the quadratic equation 

\vl (Ai) ^ + p[ (Ai) A + p 2 (Ai) = 0. (32) 

In the case at hand, the common value of Ai and A 2 lies on the unit circle. In order 
for the perturbed values of the Floquet multipliers to lie off the unit circle (and therefore 
imply instability), it is easy to verify that it is necessary and sufficient that Pi/Xi have a 
nonvanishing real part. Thus if we define 

ot = A/Ai, (33) 

we have the following criterion: 

Proposition 2 Either a is pure-imaginary and we infer stability (to leading order in e), or 
Re a 7^ and we infer instability. 

The magnitude of the real part of a is also related to the growthrate of the instability. If we 
define an instability increment 

A = |A| — 1, (34) 
then A = eRea and the growthrate is equal to A/27T, to leading order in e. 

The long calculations that lead to the coefficients appearing in equation (32) are carried 
out in the Appendices. In the notation employed there, equation (32) therefore takes the 
form 



a 2 - <J Jn + J 22 + — v + u 2 ] \ ol + 



J 11 + 2irva 1 /iJ, J 12 

J 21 J 22 + 2-Kva 2 / n 



0, (35) 



In equation (35), a has the meaning of Proposition 2 above, and the symbols J^- are defined 
in Appendix A. 2 (equation A15). There are obvious modifications of this formula if A^ = A; 
instead of Ai = A 2 . 
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3.3. The Resonant Cases 

The resonant cases for e = (circular streamlines) are those parameter values (/j,, rj) 
such that ujj — lo\ — k, where k is an integer. We'll find that these can be written in the 
form \i — f (77) (e.g., equation (37 below). Since e = 0, the \i values in question are those 
that were designated /x in equation (24) above. We no longer need the designation /i and, 
in the relations below, use the symbol ji in its place. 

If k 7^ ±2 the matrix J is diagonal (see equation A25 below) and equation (35) has the 
roots a — in + 1-aivtjj\l 1 \i and a = J22 + 27riz/u; 2 //i where <7j = iuj for j = 1, 2, 3, 4). It is easy 
to check that these diagonal entries are pure- imaginary (cf. equations A26 and A9 below) 
and therefore, in accordance with Proposition 2, there is no instability to leading order in e. 
We therefore now consider the only cases (k = ±2) that can lead to instability to this order. 

Recall that the original parameters of the problem are e - representing the departure 
from axial symmetry of the undisturbed streamlines - \x = k^/ko - representing the vertical 
wavenumber -and rj = k B. The auxiliary parameters m = fir] and q = \J,\f \ + rf are 
introduced to simplify the notation. With the frequencies taken in the order 

(ui, uj 2 , uj 3 , u 4 ) = {ji + q,-n-q,n-q,-n + q) (36) 

the replacement ix — > —\i results in the same frequencies in the opposite order. We may 
therefore assume without loss of generality that ix > 0. Further scrutiny of the formula (36) 
shows that we need only consider the following four, distinct, k — 2 resonances. 



3.3.1. Case 1. uji-uj 2 = 2 



The resonant modes are those that reduce, when 11 = 0, to the purely hydrodynamic 
modes. This case therefore represents the modification, due to the presence of the vertical 
magnetic field, of the universal elliptical instability. In this case lo\ = —002 — \x + q — 1, 
implying that 

1 



ix = . (37) 

This ratio changes from 1/2 at rj = to as i] — > 00. Evaluating the integrals defining (jij^j 
(equation A15 below) and solving equation (35) above, we find 



a 



— 7T 



2^ 

/i(l + /i) 

A 4 



(38) 
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This has a maximum instability increment (when v — /i 2 (1 + /i) /2) given (suppressing a 
factor of e) by 

®max (1 ~l~ A 4 ) • 

In the pure-hydrodynamic limit rj = we find /i — 1/2 and therefore a max = 97r/8. Since 
the growthrate is given by ae/2ir, this gives a maximum growthrate of 9e/16, in agreement 
with previous results obtained by other methods (Waleffe 1990; Goodman 1993). In the 
limit r) — > oo this maximum instability increment tends to the finite limit tt/2, about half 
its value in the pure-hydrodynamic limit. 

This instability has a bandwidth {y + — u_) e which is, for given e and rj, the length of 
the /x-interval for which the unperturbed configuration is unstable. It is determined by the 
values of v that make the real part of a vanish. These may be read off equation (38) above: 

v + = fj, (1 + fj) (1 + 3fi) /4, v_ = -n(l- fi 2 ) /4. 

In the limit 77 — > (the pure-hydrodynamic case) these give v + = 15/32 and V- = —3/32. 
These values of v± can also be inferred from Waleffe's treatment of the pure-hydrodynamic 
case. In the limit of large magnetic parameter 77, the width of the band tends to zero. 



3.3.2. Case 2. ^-^ = 2 

The resonant modes consist of a hydrodynamic mode and a purely magnetic mode (one 
frequency would vanish if r\ = 0). In this case q — 1, implying that LO\ — /i + 1, cu 3 — fj, — 1 
and 

V = -7— 2 - (39) 

v 1 + v 



Thus the ratio fi changes from 1 at rj = to as rj — > 00. This represents a 'mixed 
mode', i.e., the resonance is between a purely hydrodynamic mode and a purely magnetic 
one. Evaluating the integrals and solving the quadratic (35), we find 



a 



= mv (^j- - 1 + fi 2 ^j ± \[D 



(40) 



where 



D=-n 2 



If D < then a is pure-imaginary and the unperturbed configuration is stable, so instability 
prevails if and only if D > 0. 
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Instability can indeed occur and has its maximal increment when v = /i 3 (1 — /i 2 ) /2. 
This maximal increment is (except for a factor of e) 

(Re«)_ = ^(l-/,) 2 . 

This ranges from when 77 = 0to7r/2as?7^oo. 

The width of instability band can be calculated as in the preceding case by finding the 
values of /i for which D = 0. These are 

v + = i/i (1 - ^ 4 ) and v_ = ^ (l - fi 2 ) (3/i 2 - l) . 



3.3.3. Case 3. u 4 - lu 3 = 2 



These are the purely magnetic modes, which play no role if f] = 0. For this case we have 
q — /i — 1 implying, for fixed r] that 

" = /tA r (41) 

a/1 +77 2 - 1 

Since /x cannot exceed 1, this resonance can only occur for sufficiently large values of rj, 
namely 

T) > V3. (42) 

For a given wavelength k this would require a sufficiently large magnetic field B. When 
this condition is satisfied, we have uj± — 1 = The formula for the instability increment 
a becomes 



a 



7T 



5 (l-<0' 



7T 



h At (1 - ^) 

A 4 



(43) 



The maximum instability increment for given e is (again suppressing a factor of e) 

O^max ^ (1 AO j 

which occurs for z/ = — ii 2 (1 — //) /2. It vanishes in the limit 77 = and tends to ir/2 as 
rj — > 00. The upper and lower edges of the band of instability are expressed by the formulas 
H — /Jo + where z/± may be determined from equation (43) by requiring a to vanish. One 
finds 



X 1 

-/x (1 - n) (1 - 3/i) and i/_ = --/j (l - ii 2 ) 



The bandwidth is therefore 



v + -v_ = ^/i(l - /i) 2 . 



This vanishes both in the limit as 77 — > and in the limit as rj — > 00. The maximum 
bandwidth occurs when /j = 1/3 or 77 = \/l5. 



-14- 



3.3.4- Case 4. ui - 004, = 2 

For this fi = 1. It is clear from the expressions for the matrix J that the latter vanishes 
in this case. This leads to pure-imaginary values of a and therefore there is no instability 
associated with this resonance. 



4. Numerical results 

In this section we present a selection of numerical results. These are obtained by inte- 
grating the system (20) to obtain the fundamental matrix solution $ (t, e, fj,, if) and evaluating 
it at t — 2n to get the Floquet matrix M (e, /i, 77, ). For fixed 77, the eigenvalue of maximum 
modulus is found as a function of E (rather than e) and /i, and the regions of the E/i plane 
where this maximum modulus exceeds one are distinguished. We have carried this out to 
E — 1.6 (e = 0.4875) in the figures although there is no limitation on the size of E, or of e, 
in this method. 

In Figure 1, we have taken the magnetic parameter equal to zero in the left-hand panel, 
so this represents the purely hydrodynamic case studied originally by Bayly (1986) and 
Pierrehumbert (1986) and subsequently by others. For 77 > 0, a mixed mode of interaction 
involving both hydrodynamic and hydromagnetic modes comes into existence, and this is 
shown in the right-hand panel of Figure 1, where 77 = 1. This is too small for the remaining 
leading-order instability to appear. That remaining instability, which represents a resonant 
interaction between two modes that owe their existence to the presence of the magnetic field, 
is indicated in Figure 2 for a magnetic parameter rj = 2 (right-hand panel), slightly greater 
than the minimum value (y/3) for the existence of this instability. 

The asymptotic formulas imply that the maximal growthrate (or equivalently the max- 
imal instability increment A) for each of the wedges of instability tends to a fixed value as 
the magnetic parameter 77 increases. This is illustrated in Figure 3. However, the asymptotic 
formulas for the growthrates are less accurate than those for the stability boundaries, for the 
larger values of E. 

In identifying these tongues, we have made repeated use of the asymptotic formulas 
presented earlier. The numerical procedure also picks up some further tongues, related 
to higher-order resonances, that are excluded by the procedure leading to the asymptotic 
formulas. These we have mostly ignored on the ground that they are too weak and occupy 
too small a region of the parameter space to be significant, but we show one such resonance 
tongue for two values of rj in Figure 4 
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Fig. 1. — On the left is the case r\ = 0, i.e., the purely hydrodynamic case considered earlier by 
several other authors. It is shown here to contrast it with the case when the magnetic field is not 
zero. The lines of asterisks indicate the same stability boundaries obtained from the asymptotic 
formulas /i = fio + v±e; this approximation is seen to be quite good. It is typical of the other cases 
considered. On the right is the case r/ = 1, showing both the effect of the magnetic field on the 
hydrodynamic mode and the existence of a new mode of instability that is due to the presence of 
the magnetic field. 
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Fig. 2. — On the left is the case r] = 2. Three regions of instability occur here. The lowest of these 
is the modification of the hydrodynamic- mode instability indicated in Figure 1. The middle region 
refers to the mixed hydrodynamic-magnetic mode. The uppermost region, very thin and labeled 
MAGNETIC MODE, exists only for values of r? exceeding so is poorly developed for this value 
of 7]. On the right is the case rj = \/TE ~ 3.873, for which the width of the uppermost, purely 
magnetic, mode band is at its greatest. Note however that the vertical scale is compressed relative 
to the diagram on the left, exaggerating the widths of these bands by about a factor of two. 
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T1=0, E = 1 .3 11=3.873, E = 1.3 




0.46 0.5 0.54 [I 0.58 0.62 0.18 0.19 0.2 1^ 0.21 0.22 0.23 

Fig. 3. — The left-hand panel shows the instability increment A as a function of fi for a fixed 
value of the magnetic parameter (77 = 0) and the ellipticity (E = 1.3). The right-hand panel 
does the same except that r) = y/TE. 



11=1 r| = 3.8730 




Fig. 4. — These show higher order wedges of instability. Both are cases of resonance between 
the two hydrodynamic modes but with u>i — uj 2 = 4, rather than 2. The left-hand panel shows 
the instability wedge when 77 = 1 , the right-hand panel shows the corresponding wedge when 
7] = Vl5. 
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5. Upper bound on fieldstrength 

The results of §§3.3 and 4 show that magnetoelliptic instability occurs in vertically 
unbounded systems whatever the fieldstrength, and that the growth rates are relatively 
insensitive to the magnetic fieldstrength parameter 77. The loci of instability in the (E, //) 
plane, however, do depend on 77. As r\ — > 00, the resonant /i, in all three cases, scales as rf 1 
and the magnetic tension parameter m — > 1. This is a consequence of the resonant character 
of the instability. Since the mode frequency must be close to the rotational frequency, the 
Alfven frequency m cannot be too large. As B — > 00, k 3 must approach zero. 

In a system of finite vertical thickness H, k 3 cannot drop below x/H, where x is a factor 
of order unity. Therefore, if va exceeds VtH / x, we expect the magnetoelliptic instability to 
be suppressed. A precise upper bound on va follows from eqn. (A3) and application of the 
resonance conditions; equations (37), (39), and (41). In Cases 1 and 2, instability requires 
m 2 < 1. Case 3 requires m 2 < 3. Therefore, magnetoelliptic instability requires va < 
V?>VLH/x (there is a correction of order e due to the finite bandwidth). The corresponding 
fieldstrength is comparable to the maximum value of the field at which the magnetorotational 
instability can operate (Balbus & Hawley 1998). 

6. Discussion 

We have explored the effect of a uniform, vertical magnetic field on the stability of 
planar, incompressible flow with elliptical streamlines in an unbounded medium, in the 
approximation of ideal magnetohydrodynamics. In the absence of magnetic fields, flows 
with elliptical streamlines having ellipticity parameter e [see equation (22)] are known to 
be unstable to perturbations with wavevectors that are transverse to the plane of the flow 
(the "elliptical instability"). Our first conclusion is that this elliptical instability persists in 
the presence of a vertical magnetic field: the latter decreases the maximum growthrate but 
fails to suppress the instability, no matter how large the magnetic-field parameter becomes. 
It can be compared with the conclusion of Kerswell (1994) that a toroidal magnetic field 
has a stabilizing influence. Kerswell's analysis holds for small e only and shows that the 
growthrate decreases with magnetic field in that limit. Our result, which holds for a vertical 
magnetic field, shows a similar trend for small values of the magnetic-field parameter, but 
this trend never results in complete stabilization of the elliptical instability with increasing 
magnetic field. 

A second conclusion is that there are further instabilities associated with the presence 
of the magnetic field. One of these, for which the eigenvector is a mixture of hydrodynamic 
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and magnetic modes, occurs for all values of the magnetic field parameter. Another, for 
which the eigenvector is a combination of magnetic modes only, sets in for values of the 
magnetic-field parameter exceeding a certain threshold value (77 > V3). In all three of 
these instabilities, for large magnetic fields, the wave vector makes only a small angle with 
the plane of the unperturbed flow, reflecting the familiar tendency for dynamics to become 
nearly 2-dimensional in a strong, well ordered magnetic field. This is reflected in Figure 2, 
which shows that as rj increases, the unstable wedges are pushed to smaller /i. Although 
the unstable fraction of the (E,/i) plane decreases with increasing r\ (except for a very 
slight maximum at r\ ~ 2.18, reflecting the onset of instability between magnetic modes), 
the separation between the unstable wedges also decreases. While the nonlinear evolution 
of the unstable system is beyond the scope of this work, the destabilization of a nearly 
continuous swath of parameter space may have consequences for the interactions between 
unstable modes. In all three cases, the maximum instability increment tends to e7r/2, i.e., 
the maximum growthrate of the unstable modes tends, in dimensional units, to efi/4. 

All three of these instabilities may be relevant in accretion-disk settings. In systems of 
finite thickness H, however, the instability is suppressed if the Alfven speed va exceeds a 
critical value of order QH. Magnetorotational instabilities are quenched at approximately 
the same fieldstrength (Balbus & Hawley 1998). The growth rate of magnetoelliptical insta- 
bilities is smaller than that of magnetorotational instabilities by a factor of order e, and thus 
they are not necessarily the primary instability in magnetized disks. They may well play a 
secondary role by breaking up eddies or vortices generated by other mechanisms. 

Magnetoelliptic instabilities may also occur in the inner parts of barred galaxies, in 
which the gas flow is slightly elliptical and the magnetic field, at least in the Milky Way, has 
a vertical component (Morris and Serabyn 1996). In such settings, the instabilities could be 
a source of turbulence, possibly affecting the mass supply to a central compact object. 

We are happy to acknowledge the referee for useful comments. Material support for 
this work was provided by NSF grants AST-0098701, AST-0328821, PHY-0215581, and the 
Graduate School of the University of Wisconsin, Madison. 



Appendices 



We carry out calculations leading to the coefficients appearing in equation (32) in a series of 
steps. Since p(X) = \M — XI\, with M given by equations (25) and (26), we begin with the 
expression for M. 



A. The Expansion for M 

A.l. Zero-order Problem and Solution for M 

The matrix D (t, e, ji) has only two entries that depend on e: D n and D 2 \- If, therefore, 
we write the Taylor expansion 



D (t, e, n) = A) (t, fJ>) + eD £ (t, n) + 



(Al) 



then 



D (t, n) = D (t, 0, /x) and D e (t, fj) = D e (t, 0, //) , 
where D e = dD/de. For D we find the constant matrix 

/ -2 im \ 
2/i 2 

im 



D (t,fi)=D {n) = 



V 






im 



im 






im 





(A2) 



Its eigenvalues are 

<Ti = i (// + g) , cr 2 = -i + q) , <7 3 = i (/j, - q) , <r 4 = -i (/x - g) , 



(A3) 



where g = ^/// 2 + to 2 = /ia/1 + ^ 2 . These are distinct and nonzero as long as /i ^ and 
77 7^ 0, which we assume to be the case. The first two correspond to "hydrodynamic modes" 
since they reduce, when r\ = 0, to the eigenvalues of the purely hydrodynamic case. The 
second two refer to "magnetic modes" since they are zero in that limit. These are all of 
stable type, corresponding to frequencies Uk, k = 1, ... ,4. Regarding the matrix D t , one 
can easily work it out from the expression (21) above: all its entries vanish except (D e ) u 
and (D e ) 21 . One finds 



(A4) 



and 



(D £ ) 21 = /i 2 (1 - /i 2 ) (e 2U + e~ 21t - 2) + 2/iz/. (A5) 
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tFrom the matrices D and D e we can construct the matrices M (fx) and M e (fx) needed 
in the formula (26) for M x . For M we simply have exp 27rD . For M e we proceed as follows. 
On the finite time-interval [0, 27r] we may write 

$ (t, e, fx) = $ (t, fx) + e$! (t, fx) + . . . , $x (0, /j) = 0. 

Substituting this in the differential equation (20), expanding to first order in e, using the 
variation of constants formula (cf. Coddington and Levinson (1955)) and setting t = 2n, we 
get 

M(e,n)=M (fi)(l + e^ % x (s, n) D £ (s, n) $ (s, ds^j . (A6) 

This expresses the Floquet matrix correctly to linear order in e, and this will turn out to be 
sufficient for our purpose. The formula above identifies M e (0,/i): 

M e (0, fi) = M (fj.) [ ^o" 1 (s, fi) D e (s, fx) $ (s, fx) ds. (A7) 
Jo 

We next proceed to simplify this expression. 

A. 2. A further transformation 

The characteristic polynomial given in equation (27) above is the same in any coordinate 
system, so we shall choose one to simplify the unperturbed Floquet matrix M (fx). 

If fx ^ and rj ^ the eigenvalues {o"fc} given by equation (A3) are all distinct, so 
the eigenvectors are linearly independent and the matrix T (fx) formed from their columns 
diagonalizes D : 

D = diag(<7i, (72,(73,(74) , 

where the tilde indicates the transformed matrix: D = T~ 1 DT. We shall need to know T 
and T~ l explicitly. It is a straightforward matter to show that any eigenvector of Dq must 
have the structure (up to a constant multiple) 

\ 

- (a 2 + m 2 ) /2 
im 

im (a 2 + m 2 ) /2a J 
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Substituting the particular values of a given in equation (A3) gives the four columns of the 
matrix T, and from this we can construct its inverse. One finds 



T 



/ <7i a 2 

im 
y mix 



a 3 <t 4 \ 

-iix,o 3 iixa^ 

im im im 
-m\x mix 



-m\x J 



and 



T -i 



Amq 



f —im m/fi as icrs/pL ^ 

im mjix cr 3 —ia 3 /fi 

im —m/fi —<7i —io\jix 

\y —im —m/ji — u 1 iai/fi J 



The matrices T and T 1 depend on ix and on rj through the parameters 

m = fir] and q = a/ li 2 + m 2 = /i\A + rf. 



(A8) 



(A9) 



(A10) 



In place of equation (A6) we now obtain 

M = M (/ + ej) , 



where 



/•27T 

J(//) = / $ _1 (t,/i)D 1 (t,ii)^(t,ii) dt. 
Jo 



(All) 



(A12) 



Because the eigenvalues {at} are distinct, the matrix <f> = exp|_D ^| takes the simple, 
diagonal form 

$ (t) = diag (exp (o^t) , exp (a 2 t) , exp (<j 3 t) , exp (<7 4 i)) (A13) 
The ij entry of the matrix E>i is (since _D e has only two nonzero entries) 



{D,) t =T l3 {{T-') tl {D e ) ll + {T^) i2 {D e ) 2l ). 



(A14) 



As a result, the ij entry of the matrix J providing the leading-order perturbation of the 
Floquet matrix is 



2?r 



= (cr 3 -CTi)* 



(A) 11 (*) * 



/•27T 

+ / e^-^(D e ) 21 (t)dt. 

Jo 



(A15) 
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This enables us to find M € (0, /x). 

For the matrix M (//) we have the expression 

M (/i) = M (0, /x) = diag (Ax (/x) , . . . , A 4 (/x)) (A16) 

with Afc = 27ro"fc. Recall that in order to construct M x we need also the derivative of this 
matrix with respect to /x, 

M' (/x) = M M (0, /x) = diag (Ai (//),..., A 4 (/x)) . (A17) 

According to equations (A3) and (A10), each eigenvalue of .Do is linear in /x. Therefore 
°1 (/•*) = (aO //•*• Since \ k = exp (2nak), we have 

X' k (/x) = A fc (/x) 27rcT fc (/x) //x = \ k 27riu k /fi. (A18) 

The formulas of this section allow one to determine the matrices Mo and Mi. To 
produce from these the coefficients pj (A) appearing in equations (31) and (32) above, we 
need formulas for the derivatives of a determinant. These are presented in Section B below 
and applied in the following section. 



A. 3. The Expansion for p (A, e) 

We can now find the required expansion for p (A, e) by identifying the matrix A of Section 
B below with M — XI and the coefficients q k with the coefficients p k (A) of equation (31). We 
obtain these coefficients by writing a k = \ k — A, A' kl (0) = ( M 1 and A kl (0) = 2 ( M 2 J . 
This gives for p 1 the following expression (with n — 4): 



Pi (A) = ( 



Mi 


) n (A 2 


- A) (A 3 


- A) (A 4 


-A) 


(A19) 


Mi' 


) 22 (A. 


- A) (A 3 


- A) (A 4 


-A) + 




Mi' 


) (A, 

33 


- A) (A 2 


- A) (A 4 


-A) 




Mi' 


) (Ai 
'44 


- A) (A 2 


- A) (A 3 


-A), 





and there is a similar, lengthier expression for p 2 obtained by making the corresponding 
substitutions in equation (B5) below. 
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The development thus far has required no assumptions regarding the multipliers {A^}. 
We now suppose that Ai = A 2 . It is then clear from the expression above that pi (Ai) = 0, 
as asserted in Section 3.2. The coefficients appearing in equation (32) are now easily found 
to be 

^(Ai) = 2(A 3 -A 1 )(A 4 -A 1 ), 



and 



p'i (Ai) = - { (Mi) u + (Mi)J (A 3 - AO (A 4 - Ai) 

'Mi 



Mi 



Mi 



Mi 



P2(Ai) = 

Equation (32) therefore takes the form 



12 



22 



(A 3 - Ai) (A 4 - Ai) . 



(A20) 
(A21) 

(A22) 



Mi 



ii 



M x 



Mi Mi 

21 V / 22 



12 



0. 



(A23) 



We note that the calculation of the perturbation of Ai to first order in e, which requires 
expanding p to second order, requires the expansion of the Floquet matrix M (e) = M + 
Me H only to first order. 

We have assumed that the coincident roots are the first two, Ai = A 2 . If instead we 
should have A& = A;, equation (A23) is modified by the replacement (1,2) — > (k,l). 

Equation (A23), together with equations (All), (A16), (A17) and (A18) leads to equa- 
tion (35) of the text. What remains is to evaluate the integrals defining J, and we now turn 
to this. 



A. 4. Calculating the Elements of Mi 

By equations (26) and (A7) above (see also equation A12), the matrix Mi is given by 
the formula 

M = M (aO J + z/M M (0, /x) (A24) 

where the entries of J are given by equation (A15). From equations (A4) and (A5) it is a 
straightforward matter to carry out the integrations. We'll use for T the matrix given above 
in equation (A8). Since for this matrix = <x,- = iuij, the formula for the entries of J 
becomes 

Jij = ° 3 { (T- l ) tl ^ e^-^ (-D e ) ii (0 dt + (T- 1 )^ ^ e<"-* (A) 21 eft} . (A25) 
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For the diagonal entries the exponential factors in the integrand reduce to unity and 
one finds 

Jji = -47r/i 2 (1 - /i 2 ) a 3 (T- 1 ) . 2 , j = 1, 2, 3, 4. (A26) 

For the off-diagonal entries, the formulas may be found generally, but we need them 
only in the resonant cases where, for some pair of indices <J{ — <?j = ki for a non-zero 

integer k. 5 It is clear from the formulas (A4) (A5) and (A25) that only resonances with 
k = ±2 contribute off-diagonal terms to leading order in e since for any other choice of k the 
integrals vanish: for k ^ 2 the matrix J is diagonal. 



B. Determinantal Derivatives 



Consider an n x n matrix A (e) having the properties that it is a smooth function of e 
and is diagonal at e = 0: A (0) = diag (ai, 02, ... , a n ). 



We need the coefficients in the Taylor expansion of det A (e) = q (e) : 



Q ( e ) = <?o + qit + <?2e 2 + 



where 



9o = |A(0)|, qi 



d\A> 



de 



=0, <?2 



ld 2 \A{e) 
2 de 2 



e=Q, 



(Bl) 



:b21 



Straightforward applications of the formula for the derivative of a determinant show that 



?i = A' 11 a 2 a 3 ---a n + A' 22 aia3---a„H h A'^a^ ■ ■ • a n _i, 

• • -a n H h -4L a i°2 • • -On-i] 

^13 



g 2 = 




a 2 a 3 


+ 


A'n 
A' 

A 21 


A' 

A 1 

^22 


+ 


■ ■ + 


A 



a 3 a 4 • • • a n + 



A' A' 

■"■31 ^33 



a2<24 • • • a n 



A 1 , 



n(n— 1) 



(n— l)n 
A' 



a\Q2 ■ ■ ■ a n -2, 



(B3) 

(B4) 
(B5) 



where the terms involving two-by-two determinants represent the sum over k < I of the 
product of the {%}, and a\ omitted, with the determinant 



A' A 1 

^kk A kl 

A' A' 

^-ik A u 



and all derivatives are evaluated at the origin. 



5 Recall that at ^ <Tj for any pair i, j so k = is excluded. 
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